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We apply the density matrix renormalization group (DMRG) method to a non-equilibrium 
problem: the asymmetric exclusion process in one dimension. We study the stationary state of 
the process to calculate the particle density profile (one-point function). We show that, even 
with a small number of retained bases, the DMRG calculation is in excellent agreement with 
the exact solution obtained by the matrix-product-ansatz approach. 
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§1. Introduction 

Many phenomena in the vast realm of non-equilibrium physics can be viewed as the asymmetric 
(simple) exclusion process (ASEP for short). Models of the ASEP in one spatial dimension 
have been extensively studied since a remarkable papeii^ appeared, which proposed the so-called 
matrix product ansatz (MPA). For a class of models, the exact stationary probability distributions 
have been obtained by the MPA. However, there are many models whose stationary probability 
distributions cannot be obtained by using the MPA. Alternative approaches are required to study 
general models of the ASEP. 

There are two types of models for the ASEP: discrete-time models and the continuous-time mod- 
els. The discrete-time models which we treat in this paper are more general than the continuous- 
time ones because the latter are realized by taking suitable limits of the parameters in the former. 

The density-matrix renormalization group (DMRG) method invented by S. R. White is an ap- 
proximate numerical scheme for diagonalization of Hamiltonians, which has been applied to various 
one-dimensional quantum systems, with results demonstrating its surprising efficiency As pointed 
out by T. Nishino, the DMRG is also applicable to two-dimensional lattice statistical models for 

* e-mail: hi6ida@godzilla.phys.sci.osaka-u.ac.jp 
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diagonalization of the transfer matrices.^) Since the time-evolution operator in the discrete-time 
ASEP can be regarded as a transfer matrix, we can naturahy expect apphcabihty of the DMRG also 
to the ASEP. In the actual implementation, however, we should consider the following two points. 
First, we should treat a non-symmetric transfer matrix. Second, the role of the wave function ^ is 
different; in the ASEP not \ip\'^ but ip itself is the probability distribution. 

In this paper, by properly incorporating the above two points, we apply the DMRG to a model 
of the ASEP to show that the DMRG is an efficient numerical approach also to non-equilibrium 
problems. 

§2. Model 

Let us explain the model which we treat. It is the model of stochastic process of hopping particles 
on a chain of size A'^. For technical simplicity, we will assume to be 4 x (integer) -|- 2. The time 
is taken to be discrete. Each site can have two states: either empty or occupied. Particles hop to 
the right-nearest site with probability p and to the left-nearest site with probability q. The chain 
are coupled to particle reservoirs at the both ends: if the leftmost site is empty (resp. occupied), a 
particle is added (resp. removed) with probability a (resp. 7); if the rightmost site is empty (resp. 
occupied), a particle is added (resp. removed) with probability S (resp. /3). 

Time evolution of the model consists of the following two steps. 

[A-1] 

For all odd-even pairs of sites {{2j — 1, 2j)} (j = 1, 2, 3, A^/2 — 1) their states are updated simul- 
taneously. States of the leftmost and the rightmost sites are updated in this step. 



[A-2] 

For all even-odd pairs of sites {(2j, 2j + 1)} {j 
taneously. 



= 0, 1, 2, A^/2 — 1) their states are updated simul- 



We introduce a transfer matrix T which describes the discrete-time evolution of the probability 
distribution. Corresponding to the two steps of the time evolution of the model, T consists of two 
transfer matrices Ti and T2: 

T = TaTi , (1) 

where Ti and T2 correspond to the first step [A-1] and second step [A-2] of the time evolution, 
respectively. Explicitly, 

Ti = £ (g) Ti,2 ® T3,4 • • • T^-3,iV-2 7^ 
T2 = To^i 72,3 ■ • • '?iV-2,Ar-l , 
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where £ and TZ are the left boundary transfer matrix and the right boundary transfer matrix, 
respectively and is a local transfer matrix acting on the pair of the sites i and i + 1. 

In order to represent explicit matrix elements, we introduce a state variable in at each site n as: 



i„ = (empty); 
Then, the local transfer matrix T is given by 



1 (occupied). 



(3) 
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(4) 



where the element {Tn^n+i)(^f^ j (j^ j represents the probability for the transition from the 
state {in-,in+i) to the state {fnifn+i}- The boundary transfer matrices are 



C 





/ 1 - a 



1 



a 



7 

1-7 



,7^ 




(5) 



The stationary state of the model has three phases: the high density phase, the low density 
phase and the maximal current phase. They are separated by the boundary-induced phase transi- 
tionsJ3).lll)|) 

§3. Method 

In this article, we concentrate on the density profile (one-point function) (0 < n < — 1) 
in the stationary state, which can be calculated in two steps as follows. 



[B-1] 

Obtain the right largest-eigenvalue eigenvector IP) of T which satisfy 

T\P) = \\P), 



(6) 



with A(= 1) being the largest eigenvalue of T. This \P) is the stationary state whose wavefunction 

(io, ^1, ^2, ■ ■ ■ , ^n, ■ ■ • , iN-i \P) gives the stationary probability distribution. 

[B-2] 

Calculate {in): 

,. > _ E{j}^n(^0,il,'i2,--- ,^n,--- ■,iN-2,iN-l\P) 



E{i} («0,«1,^2, 



,iN-2-,iN-l\P) 



(7) 



where Yl,{i} denotes summation over all states {ig, ii, ■ • • , iAr_2, ^Af-i} of sites. 
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We should note that to calculate the quantity (i„), we do not use the left largest-eigenvalue 
eigenvector of T. This contrasts with the case of the classical statistical system, where the one- 
point function is expressed as {Q\in\P) / {Q\P) with {Q\ being the left largest-eigenvalue eigenvector 
of T and in the local density operator. 

Let us explain how we can incorporate the DMRG in these processes [B-1] and [B-2]. The 
method we employ is essentially the finite-system-algorithm versioni^ of the "classical" DMRG.B* 
However, since our transfer matrix T is not symmetric, we need a non-trivial modification in the 
actual application of the algorithm: Unlike the ordinary DMRG, we should deal with "asymmetric 
density matrix" 

To clarify our explanation, we treat a chain of size = 14. Suppose that, at some iteration 
stage in the DMRG, the renormalized transfer matrix T is expressed as in Fig. 1. For such form 
of T, we adopt the notation T^^^xT^, where T4 is the renormalized left block transfer matrix of 
4 sites, "x" represents the local transfer matrix and Tg is the renormalized right block transfer 
matrix of 8 sites. Just as in the original finite-size algorithm of the DMRG, we consider the set 
{TIx^ xT^-k~2}k=2,4,---,N-4 to refine iteratively the largest-eigenvalue eigenstate. 

Hereafter we use Greek indices to represent the block state. By the power method (i.e., multi- 
plication of T sufficiently many times), we get the right eigenfunction ipifj-s; ^3, «4, ^5, ie! ^^7) and the 
left eigenfunction (t>{^J■3',i3,^4,i5,^6',l^7) belonging to the largest eigenvalue of T = T^x^xTg. We 
then form the "density matrix"&'&'0)'0)'lli) for the left 5 sites p''^ = {p''^{p3,is,i4; fJ-'^ji'^ji'^)} as 

/'5(/^3,«3,^4;^3!^3>^4) = i^il^S'^ ^3, H, k, k] '^7)HfJ'3i ^5, k] l^v) ■ (8) 

is, H, 1^7 

Since (f> tlj/m general, for asymmetric T, the matrix pi^ is also asymmetric: pi,{p3, ^3,^4; /^s; ^3) ^4) 7^ 
/95(/i3, ^4; ^3,13,14). As in the ordinary DMRG, we must perform eigenvalue decomposition of P5 
to choose a truncated basis set for a new 5-site block. In our case, we should obtain both of the 
right and left eigenvectors for this purpose. Using the obtained transformation matrices for the 
new block states, we can renormalize the left block transfer matrix of 6 sites T^x^ to get Tg. All 
the necessary block-making procedure can be done in a similar fashion. 

By {T|^, ^-^ ^-^ ^ } we denote the components of the pj-th right eigenvector of pj (density 
matrix for left j sites) with eigenvalue w^^. : 

A'j-2,«j-2,«j-l 

where we set pi = iq. For the density matrix for right j sites, we denote its right eigenvector 
with eigenvalue w'l. by W , whose components {R'l^i^ . +1 u 2) ^ satisfy 

E P'j(^N~j,i'N~j+l,''j~2;iN~j,iN~j+l,l^j-2)Rli^_.^i^_^^^^^^_^)^^^ 

«JV-j.*JV-j + li^'i-2 
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where we set ui = in-i- We can then express the stationary-state vector \P) in eq.(^) as 

(^o,n, • • ■ ,k2,hz\P) 

53 



?7 



/?' /?^ /?'^ 



The one-point function can be obtained as follows. For concreteness, we explain calculation 
of (^4) and (is). 



[C-1] 

Perform partial summation 

X! *' '^\^o,k, ■ ■ ■ ,iN-2,iN-i\P) = g{u,ib) , 

ii} 

where J2{i}^^ means that we do not perform the summation over Z4,i5. 



(12) 



[C-21 



From the definition of (in) (eq.(|^) 





^g{u = I, id) 


/ 














XI 5(^4, ^5 = 1) 


/ 






14 







(13) 
(14) 



we have 



By introducing the following quantities, 

iii,n2,«i3 

E 49,no,.3),.5^3(-3)-^5(-5) 



19,110, i^3 



(15) 
(16) 
(17) 
(18) 



5'(h,«5)= E Sl{^^)i:{^i'i]i'i,iA,i^,iQ]Vj)STi{vj). (19) 

Eq.(p!9[) means that for calculation of (^4) and (is) we do not need the full probability distribution; 
we need only 53(^3), ^/'(^a; ^3, ^4, ^5, ^e; ^i): Sj{v'j). Accordingly, we need not retain full set of arrays 



{/ij_2,«j_2,ji. 



N } and {R-'f ^ | in eq.(ll) at each iteration. 



5 



§4. Results 

We present the calculated results for two cases: (A) p = 0.75, q = 0.25, a = 0.9, /3 = 0.8,7 = 
0.01,5 = 0.01, and (B) p = 0.75, q = 0.25, a = 0.5, /3 = 0.5,7 = 0.4,5 = 0.9. The case (A) 
corresponds to the maximal current phase, and (B) to the high density phase. 

In Fig. 2, we compare the DMRG calculation for case (A) with = 10 to the exact results (both 
by the exact diagonalization and by the programB using the MPA), varying the number m of the 
retained bases. It should be noted that a very small number (m = 2) of retained bases gives a quite 
accurate density profile, implying that the "weight" of the density matrix is "localized" within very 
small number of bases. The similar behaviors are seen in the DMRG calculation for N = 102 in 
both cases (A) and (B), which are given in Fig. 3 and Fig. 4, respectively. The DMRG calculations 
with a very small number (m = 2) of retained bases are also in excellent agreement with the MPA 
results. 

It has been known that for a set of parameters, the algebra associated with the MPA admits 
finite-dimensional representations.@'''lll* In the n-dimensional case, the DMRG with n retained bases 
should be exact. The parameters in the present two cases (A) and (B) do not correspond to such 
finite-dimensional cases, which means that infinite m in the DMRG is required for exact calculation. 
Our results show that even for such "infinite-dimensional cases", DMRG calculation with a very 
small number of m serves as an excellent approximation. 

§5. Summary 

In this paper, we have applied the density matrix renormalization group (DMRG) to asymmetric 
(simple) exclusion process (ASEP) which is a non-equilibrium many-body problem of hopping 
particles in one dimension. We have calculated the density profiles (one-point function) in the 
stationary states, which are in excellent agreement with exact solution. We have found that only a 
very small number of retained bases in the DMRG is necessary to attain the close-to-exact results. 
We have thus verified that the DMRG is an efficient method to study the ASEP. 

We should comment on possible extensions of the present work. 

1. Although we have considered only the one-point function in this paper, we can easily calculate 
2-point function (imin), and, in principle, general n-point function in the stationary state. This 
is because, DMRG allows us to obtain the full distribution function in the stationary state. 

2. The ASEP can be generalized in many different ways. For these generalized models stationary 
solutions cannot always be obtained by the matrix-product-ansatz approach. The DMRG 
approach used in the present paper can be easily generalized to these models. For instance, 

• Models including the processes such as coagulation, decoagulation, (pair-)creation, and (pair- 
)annihilation (We can incorporate this generalization in eq.(|^.) 

• Models of site-dependent parameters of probability p{n),q{n), where n represents a site in 
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the system 

• Models with iV-species particles (iV > 2). 

• Models of hopping particles on a ring with a defect. (The DMRG can deal with 1-dimensional 
system with periodic boundary condition.) 

Some of the above are now undertaken, whose results will be given in future publications. 
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Fig. 1. A form of the transfer matrix T = Tix'^xTg = 14)- According to Nishino's diagramja^ Q, □ and O 
represent a single site, a 3-site block and a 7-site block, respectively. 



Fig. 2. The density profile of the model in the maximal current phase (A'^ = 10, p = 0.75, g = 0.25, a = 0.9, /3 = 
0.8,7 = 0.01,5 = 0.01) 



Fig. 3. The density profile of the model in the maximal current phase [N — 102, p = 0.75,5 = 0.25, a — 0.9, /3 = 
0.8,7 = 0.01,5 = 0.01) 



Fig. 4. The density profile of the model in the high density phase {N = 102, p = 0.75, q = 0.25, a = 0.5, (3 = 0.5, 7 — 
0.4,5 = 0.9) 
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